Hydration of an apolar solute in a two-dimensional waterlike lattice fluid 
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In a previous work, we investigated a two-dimensional lattice-fluid model, displaying some wa- 
terlike thermodynamic anomalies. The model, defined on a triangular lattice, is now extended to 
aqueous solutions with apolar species. Water molecules are of the "Mercedes Benz" type, i.e., they 
possess a Da (equilateral triangle) symmetry, with three equivalent bonding arms. Bond formation 
depends both on orientation and local density. The insertion of inert molecules displays typical 
signatures of hydrophobic hydration: large positive transfer free energy, large negative transfer en- 
tropy (at low temperature) , strong temperature dependence of the transfer enthalpy and entropy, 
i.e., large (positive) transfer heat capacity. Model properties are derived by a generalized first order 
approximation on a triangle cluster. 

PACS numbers: 05.50.-|-q 61.20.-p 65.20.-l-w 82.60.Lf 



I. INTRODUCTION 

The term hydrophobicity Q refers to peculiar ther- 
modynamic properties of the transfer process of an ap- 
olar solute into water. In such a process, one generally 
observes large positive transfer free energy, large nega- 
tive transfer entropy (at room temperature), and strong 
temperature dependence of the transfer enthalpy and en- 
tropy, i.e, large (positive) transfer heat capacity. Such 
anomalous behavior seems to play a central role for rele- 
vant phenomena taking place in aqueous solutions, such 
as folding of macromolecules anclproteins, and forma- 
tion of micelles and membranes For example, bio- 
polymcrs such as proteins contain a significant fraction 
of apolar chemical groups, and it is well established that 
the efi^ective attraction between apolar groups, due to hy- 
drophobicity, gives an important contribution both to the 
folding process and to stabilization of the folded protein. 

Despite several decades of research, the theory of the 
hydrophobic effect is still incomplete. Different theoret- 
ical models have been proposed to explain the anoma- 
lous behavior of water itself and hydrophobic interac- 
tions, which have been recognized to be closely related. 
A straightforward approach relies on the application of 
molecular dynamics or Monte Carlo simulation methods 
to models with more or less realistic three-dimensional 
water geometry 0jSIE@;0- This approach is powerful, 
but has some limitations. First of all, large computa- 
tional effort is needed, and properties involving multiple 
derivatives of the free energy (such as transfer heat ca- 
pacity, whose behavior is a fingerprint of hydrophobicity) 
cannot be determined easily. Moreover, the amount of 
different interactions and geometric parameters included 
in the model makes it difficult to find out relevant phys- 
ical mechanisms which determine observable properties. 

A complementary approach involves investigation of 
simplified models IsHiliflj with fewer parameters, ge- 
ometric details, and degrees of freedom. Such models 
should better allow to trace connections between mi- 
croscopic structure and observed properties, while the 



latter can be usually analyzed in full detail, and in a 
large range of thermodynamic conditions, with relatively 
small computational effort. One of these attempts is 
based on the application of scaled-particle theory 
to hydrophobic hydration 0, These studies sug- 

gest that the hydrophobic effect results mostly from the 
small size of water molecules, and not from water struc- 
turing by the solute. A recent and interesting descen- 
dant of scaled-particle theory is the information theory 
approach by Pratt and coworkers [T^ ITU , based on pre- 
vious knowledge of water properties, such as the pair 
correlation function, obtainable either by experiments or 
by simulations. The latter approach suggests that water 
structuring induced by the solute, though existing, may 
be scarcely relevant for a description of the hydrophobic 
effect. The simplified molecular thermodynamic theory 
of Ref. is essentially in agreement with this conclu- 
sion. On the contrary, different theories stress that the 
large positive heat capacity variation, observed upon in- 
sertion of apolar solutes into water, can only arise from a 
cooperative phenomenon, that is from induced ordering 
of water molecules, so that a theory of the hydropho- 
bic effect should be based on a description of this phe- 
nomenon. It can be observed that, at room tempera- 
ture, hydration of apolar solutes is energetically favor- 
able but sufficiently unfavorable entropically, with a re- 
sulting increased free energy. A simplified way to re- 
produce this effect is for instance the one-dimensional 
lattice model by Kolomeisky and Widom extended 
also to two and three dimensions that case, en- 

tropy penalty is achieved by lowering the large number of 
possible orientations of a water molecule, if it has to ac- 
commodate a neighboring hydrophobic solute. Another 
possibility is to give water molecules a geometric struc- 
ture, as in the two-dimensional Mercedes Benz (MB) 
model, first introduced by Ben-Naim in 1971 
the latter model, water molecules possess three equiva- 
lent bonding arms arranged as in the Mercedes logo, and 
a geometric constraint (arm alignment) is required for 
bond formation. An MB model has been investigated 



2 



quite recently by Dill and coworkers by means of sev- 
eral different methods: Constant pressure Monte Carlo 
simulations 0, Eoll entr opy expansion , and integral 
equation theory |22. 01 • A coherent picture of the 
hydrophobic effect phenomenology has been worked out, 
allowing to obtain a microscopic view of several anoma- 
lous properties of water both as a pure substance and as 
a solvent. The latter studies follow the previously men- 
tioned idea of simplified models, that nevertheless are 
based on well defined microscopic interactions, that is, 
on an energy function, without previous knowledge of 
water properties. One important reason to do so is the 
need of modelling water in a computationally convenient 
way, in order to investigations on complex systems such 
as biomolecules, for which water plays a key role. 

According to the same idea, we have recently inves- 
tigated [2^ a lattice-fluid model of MB type on a (two- 
dimensional) triangular lattice, with the aim of reproduc- 
ing qualitatively the thermodynamic anomalies of pure 
water. Of course, the lattice allows important simplifica- 
tions, so that a sufficiently accurate modelling of water 
on a lattice may be quite an interesting issue. Our model 
hamiltonian includes Van der Waals interaction and hy- 
drogen bonding, when two nearest neighbor MB water 
molecules point an arm to each other. Bonds can be 
weakened by the presence of a third competing molecule, 
close to the formed bond, to mimic the fact that, if water 
molecules arc too close to one another, hydrogen bonds 
may be perturbed or broken. As far as bonding prop- 
erties are concerned, the model is equivalent to an early 
model proposed by Bell and Lavis [2g, but the weaken- 
ing criterion is different, that is, equivalent to the one 
employed in quite a recent investigation by Patrykiejew 
and coworkers |23, 13 ■ Nevertheless, in the latter anal- 
ysis the possibility of nonbonding orientations for water 
molecules is neglected, and resemblance with real water 
behavior is poor. In Ref. we observed that, intro- 
ducing such additional orientations, which account for 
directionality of hydrogen bonds, the lattice model de- 
scribes several anomalous properties of pure water in a 
qualitatively correct way. Here we extend the model to 
the case of an aqueous solution, working out solvation 
thermodynamics for an inert (apolar) solute. The anal- 
ysis is also extended to transfer properties of water in 
its own pure liquid. Our purpose is to verify whether 
this simple model, which nevertheless accounts for the 
competition between Van der Waals interactions and hy- 
drogen bonding in almost the same way as the off-lattice 
MB model, may be able to reproduce at least the main 
features of hydrophobicity. In particular, we also inves- 
tigate how the model describes the solvation process at 
a microscopic level, by comparing the average number 
of hydrogen bonds formed by water molecules (hydro- 
gen bond coordination) in the bulk or in the hydration 
shell. We shall carry out the analysis by means of a gen- 
eralized first-order approximation on a triangle cluster, 
which has been verified to be quite accurate for the pure 
water model psj . 



The paper is organized as follows. In Sec. II we define 
the model and recall the first-order approximation. In 
Sec III we introduce the thermodynamic quantities used 
to characterize the solvation process, and compute them 
in the framework of the model. In particular, we ana- 
lyze transfer quantities for an inert molecule, comparing 
them to the case in which hydrogen bonding interaction 
is "turned off" . Model predictions about solvation of wa- 
ter in its own pure liquid are also reported. In Sec. IV 
we investigate hydrogen bond coordination, drawing a 
comparison with the results of the off-lattice MB model. 
Sec. V contains some concluding remarks. 



II. THE MODEL AND THE FIRST ORDER 
APPROXIMATION 

The model is defined on a two dimensional triangu- 
lar lattice. A lattice site can be empty or occupied by 
molecules of two different chemical species, water (w) or 
solute (s). A water molecule has three equivalent bond- 
ing arms, separated by 27r/3 angles, whereas a solute 
molecule is assumed to have no internal structure. Two 
nearest-neighbor (NN) molecules of species x, y (with 
X, y = w, s) interact with an attractive energy — Cxy < 0, 
representing ordinary Van der Waals forces. Moreover, if 
two arms of two NN water molecules point to each other, 
an energy term —77 < is added, to mimic the formation 
of a hydrogen (H) bond. Due to the lattice structure, a 
water molecule can form 3 bonds at most, and has only 2 
bonding orientations, when its arms are aligned with the 
lattice. We also assume that w nonbonding configura- 
tions exist, where w is a model parameter, related to the 
H bond breaking entropy. The H bond is weakened by an 
energy term Cx?7/2 (cx £ [0, 1]) when a third molecule of 
X species is on a site near a formed bond. In the triangu- 
lar lattice there are 2 such weakening sites per bond, so 
that, when both are occupied by x molecules, the H bond 
contributes an energy —(1 — Cx)?;. The resulting water- 
solute interaction has two components: The nonorienta- 
tional Van der Waals term — e^s and the weakening term 
Cs?7/2. The latter, which is an effective 3 body interac- 
tion, can be viewed as a perturbation effect of the solute 
on a H bond between two water molecules. 

Let us write the model hamiltonian. In order to in- 
troduce the first order approximation, it is convenient to 
express it as a sum over elementary triangles: 

n = l ^Mwv., (1) 

where Tiijk is a contribution which will be referred to as 
triangle hamiltonian, and ir,ir',ir" label site configura- 
tions for the 3 vertices r, r', r", respectively. Possible site 
configurations are (see Tab.Q}: "empty" (z = 0), "bond- 
ing water" (site occupied by a water molecule in one of 
the 2 orientations which can form bonds: i = 1, 2), "non- 
bonding water" (site occupied by a water molecule in one 
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TABLE I: Possible site configurations, with corresponding 
labels (i), multiplicities (wi), and occupation numbers for wa- 
ter {n„^i) and solute {ns,i). 



config. empty 




i 1 2 3 4 

Wi 1 1 1 TO 1 

n„,i 1110 

ns,i 1 



of the w orientations which cannot form bonds: i = 3), 
"solute" (site occupied by a solute molecule: « = 4). Even 
if all configurations have unit multiplicity, except i = 3, 
it is convenient to introduce a generic multiplicity param- 
eter Wi, defined as in Tab. The triangle hamiltonian 
can be written as 

Hijk = Hijk + Hjki + Hkij, (2) 

where 

Hi;jk = -exyn.x,i"yj - vhij (1 - Cxrix.fc) , (3) 

nx,i is an occupation variable for the x species, defined 
as in Tab. while hij = 1 if the pair configuration 
forms a H bond, and hij — otherwise. Con- 
ventionally, repeated x and y indices are summed over 
their possible values w, s. Let us notice that triangle ver- 
tices are set on three triangular sublattices, say A, B, C, 
and i,j,k are assumed to denote configurations of sites 
placed on A, B,C sublattices respectively. Assuming also 
that A, B, C are ordered counterclockwise on up-pointing 
triangles (whence clockwise on down-pointing triangles), 
we can define hij = 1 if i = 1 and j — 2 and hij = 
otherwise. Let us also notice that both Van der Waals 
(— exy^T-x.^ny j) and H bond energies {—rjhij), that arc 2- 
body terms, are split between two triangles, whence the 
1/2 prefactor in Eq. On the contrary the 3-body 

weakening terms {rjhijCxn^^k/2) are associated each one 
to a given triangle, and the 1/2 factor is absorbed in 
the prefactor. Let us denote the triangle configuration 
probability by Pijk-, and assume that the probability dis- 
tribution is equal for every triangle (no distinction be- 
tween up- or down-pointing triangles). Taking into ac- 
count that there are 2 triangles per site, we can write the 
following expression for the internal energy per site of an 
infinite lattice 

4 4 4 

u ^^^^WiWjWkPijkT-Cijk- (4) 

i=0 3=0 k=0 

The multiplicity for the triangle configuration {i,j,k) is 
given by WiWjWk, where Wi is the previously mentioned 
multiplicity parameter. 

The details of the finite temperature analysis of the 
model, by a generalized first order approximation on a 



triangle cluster, follow the previous paper |25|. Not being 
interested in symmetry broken phases (ice) , we introduce 
a homogeneity condition since the beginning, after which 
generalization is trivial. The grand-canonical free energy 
per site uj = u — n^Px — Ts {^-^ and px being respectively 
the chemical potential and the density, or site occupation 
probability, for the x species, T and s being respectively 
the temperature and the entropy per site), can be written 
as a functional in the triangle probability distribution as 

4 4 4 

w = ^ ^ ^ w^WjWkPijk (iiijk +Th-ipijk^ 

i=0j=0k=0 

- 2T^WiPi \npi, 

i=0 

where T is expressed in energy units (whence entropy in 
natural units) . It is noteworthy that the only variational 
parameter in cu is the triangle probability distribution, 
that is the 125 variables {pijk}-, because pi (the site prob- 
ability) can be expressed as a marginal. According to the 
homogeneity hypothesis, we can use the symmetrized ex- 
pression 

j=0 fc=0 

The latter turns out to be convenient, in that it gives 
rise to iterative (fixed point) equations which preserve 
homogeneity. The energy part in Eq. (jSj) is given by the 
modified hamiltonian 'Hijk, defined as 

i^ijk — rlijk — Mx ^ • 

The minimization of to with respect to {J>ijk\^ with the 
normalization constraint, yields the equations 

P,:,fc=r'e-^--/^(p,;p,Pfc)'/', (8) 

where ^ is a normalization constant. Such equations can 
be solved numerically by simple iteration (natural iter- 
ation method H^). As a result, we obtain the triangle 
probability values at equilibrium, from which one can 
compute the thermal average of every observable. The 
substitution into Eqs. and ISJ gives respectively the 
equilibrium internal energy and free energy, and, by the 
way, it is possible to show that 

cj = -Tln^. (9) 

III. SOLVATION THERMODYNAMICS 

Let us first introduce the thermodynamic quantities, 
generally used to describe solvation from a macroscopic 
point of view, which we shall evaluate for our model in 
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the following. The physical process we are interested in is 
the transfer of a component (x = w, s) into a water-solute 
mixture, with solute density tending to zero. According 
to the Ben-Naim standard i this process is well char- 
acterized by variation of the pseudo-chemical potential 
fi*, related to the ordinary chemical potential jUx of the 
given component by 



+Tlogpx- 



(10) 



Let us notice that, in the latter term on the right hand 
side, the momentum partition function is missing, due 
to the fact that we are dealing with a lattice model [sOl- 
The solvation free energy per molecule Ag* is defined as 
the difference between pseudo-chemical potential values 
of a molecule x in the ideal gas phase (g) and in the liquid 
phase (/). For practical purposes, the differences between 
ideal and real gas can be generally neglected [s^l- For a 
molecule of species x we then have 



(I) 



(9) 



(11) 



where fix and fix are pseudo-chemical potentials of 
X in the liquid and gas phases, respectively. Now, if liq- 
uid and gas phase coexist in equilibrium, as usual in ex- 
periments, the ordinary chemical potentials of x in both 
phases must be equal, and by simple algebra we obtain 



As: 



-Tin 



Px 



(12) 



where px ^ and px^^ are respectively the x species densi- 
ties in the liquid and in the gas. Derived quantities, of 
interest in experiments, are the solvation entropy 



As! 



the solvation enthalpy 



dAg; 



dT 



Ah:^Ag:+TAs:, 
and the solvation heat-capacity 



Acp! = 



dAht 



dT 



(13) 



(14) 



(15) 



Let us notice that, in principle, we should distinguish be- 
tween a constant pressure derivative (as stated by defini- 
tion) and a temperature derivative taken along the liquid- 
vapor equilibrium curve. In particular, we could not even 
use Eq. (|12|l . because we would move out of the equilib- 
rium curve, at which the ordinary chemical potentials 
are equal. Nevertheless, numerical results show that 
the difference is negligible, definitely from a qualitative 
point of view. Therefore, we shall always use the latter 
temperature derivative definition in our calculations. 

Let us start studying solvation of an inert molecule in 
water. We set water parameters as follows: 77/eww — 3, 



c„ = 0.8, and w = 20. As shown in our previous analy- 
sis [23, this choice corresponds to a waterlike behavior, 
with a liquid- vapor critical point at r/e„„ « 1.18, and 
with a temperature of maximum density for the liquid 
around T/e^w « 0.67 at low pressure. Solvation ther- 
modynamics concepts introduced above are independent 
of density of components in the mixture. We choose to 
let solute density assume very low values with respect 
to water density (dilute solution limit), in order to com- 
pare with experiments. The "inert" nature of the so- 
lute is described, in the model framework, by setting to 
zero the solute-solute (ess) and water-solute (ews) Van der 
Waals interaction energies. At the moment, we also set 
the weakening parameter Cs to zero, assuming that the 
ideally inert solute does not weaken H bonds in its neigh- 
borhood. The effect of nonzero values for this parameter, 
which may describe for instance a volume interaction, 
resulting in a perturbation of the geometry (and hence 
of the energy) of H bonds, will be taken into account 
later. The temperature trends of the free energy, en- 
tropy, and enthalpy of transfer are given in Fig.QJa); the 
transfer heat capacity in Fig. ^c). In order to compare 
with experimental data [30| . all quantities are evaluated 
at liquid-vapor coexistence. For the perfectly inert so- 
lute, we have verified that concentration does not affect 
the results at all. To represent roughly the experimental 
temperature range (between 0° C and 300° C) we report 
model results between T/e„y^ = 0.65 (just below the tem- 
perature of maximum density for pure liquid water) and 
T/eww = 0.90 (about half way between the previous tem- 
perature and the critical temperature). Remarkably, it 
turns out that the model, despite its simplicity, displays 
the defining feature of the hydrophobic effect: the sol- 
vation free energy is positive and large, while the solva- 
tion entropy is negative at low temperatures and becomes 
positive upon increasing temperature; the heat capac- 
ity is positive and large, and also the decreasing trend 
with temperature is essentially reproduced. The increas- 
ing trend at high temperature is related to the the fact 
that we are approaching the liquid- vapor critical point. 
Negative solvation entropy at low (room) temperature is 
a clear indication that solute insertion into the mixture 
orders the system. The corresponding positive (unfavor- 
able) contribution to free energy compensates a negative 
(favorable) enthalpic contribution, giving rise to a pos- 
itive solvation free energy. At higher temperature, en- 
thalpic and entropic contributions change sign, but they 
still have the same trend of compensating each other. As 
observed in experiments, the model predicts two differ- 
ent temperatures Th and Ts at which the transfer en- 
thalpy and entropy vanish, respectively [see Fig. Q^b)] . 
This behavior is to be ascribed to the thermodynamics 
of H bonding and, in order to rationalize this fact in the 
model framework, let us also analyze transfer quantities 
for the case ry = 0, i.e., turning off H bond interactions 
[see Fig. ^a),(c)]. As expected, the results are similar 
in the high temperature regime, where there is a high 
probability that hydrogen bonds are broken by thermal 
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FIG. 1: (a) Solvation energies (-E/cmv) vs temperature (T/e^) for an ideal inert molecule in water at liquid-vapor coexistence: 
E = Ag* (solid line), E = T As* (dashed line), and E — Ah* (dashed-dotted line), (b) Corresponding experimental data for 
transfer of argon into water |30ll . (c) Solvation heat capacity (Acp*) vs tem per ature (r/e„„) for an ideal inert molecule in 
water, (d) Corresponding experimental data for transfer of argon into water 13(1 . Thin lines in (a) and (c) denote solvation 
energies in nonbonding water {rj — 0). 



fluctuations, whereas they change more and more dra- 
matically upon decreasing temperature, and in particu- 
lar the regions of negative transfer entropy and enthalpy 
completely disappear. This facts confirm that H bond- 
ing is the key element for system ordering, upon insertion 
of an inert molecule. Accordingly, also the divergentlike 
trend of the heat capacity upon decreasing temperature 
(related to the fact that the liquid phase is approach- 
ing a stability limit p5l |) is suppressed. The process is 
now completely dominated by the transfer enthalpy, with 
a large and positive transfer free energy, and a positive 
transfer entropy. The transfer quantities behave qualita- 
tively as those observed in solvation experiments of noble 
gas molecules in ordinary liquids pMl ISlj , and are rela- 
tively independent of temperature. In fact, with 77 = 0, 
a water molecule can be viewed as a nonpolar molecule 
with Van der Waals interaction energy e^w 

Let us now consider also the solvation of water in its 
own pure liquid. The corresponding transfer energies ob- 
tained by the model are displayed in Fig.|2Ia), where we 
have reduced the temperature interval, in order to com- 
pare with available experimental results |30| , reported in 
Fig. |2tb). In contrast to the inert molecule case, here 



the absolute values of solvation free energy and entropy 
are considerably lower and the enthalpy, rather than the 
entropy, dominates the solvation process, and all quan- 
tities are relatively independent of temperature. These 
features characterize a regular transfer process, like the 
solvation of an ordinary fluid molecule from a gas phase 
into its pure liquid phase. In this case, upon removing 
H bond interactions [thin lines in Fig. |2Ia)], very little 
changes arc observed in the solvation energies, except at 
very low temperature, where we are approaching the sta- 
bility limit for liquid water . Let us discuss two issues 
about these results. First, the fact that so little changes 
are caused by turning on or off H bonds can be ratio- 
nalized on the basis of the microscopic model interac- 
tions. The insertion of a water molecule into pure liquid 
water should imply in principle the formation of new H 
bonds, but the model is such that insertion of a new wa- 
ter molecule also weakens other H bonds in its neighbor- 
hood, and the two effects nearly compensate each other. 
Second, let us notice that solvation enthalpy decreases 
upon increasing temperature, that is, the solvation heat 
capacity is negative, in contrast with experiments. We do 
not have an explanation for this fact, but we can observe 
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FIG. 2: (a) Solvation energies (iJ/e„„) vs temperature (r/e,™) for a water molecule into pure liquid water at liquid- vapor 
coexistence: E = Ag^ (solid line), E = TAs* (dashed line), and E = Ah* (dashed-dotted line), (b) Corresponding experimental 
data [soil . Thin lines in (a) denote solvation energies for nonbonding water (77 = 0). 



that an analogous effect is observed when the model is 
reduced to describe a regular solvation process, that is 
when H bonds are turned off. This suggest that there 
is probably a limitation of the lattice description, that 
anyway has nothing to do with the peculiarities of water. 
The effect is quantitatively small, so that it is hidden by 
other large cnthalpic and cntropic effects observed in the 
case of hydrophobic solvation. 



IV. HYDROGEN BOND COORDINATION 

So far, we have always considered an ideal, perfectly 
inert solute, setting to zero all interactions with water 
(ews, Cs) and with itself (£55). Now we investigate the role 
of the Cs parameter, which represents a weakening of H 
bonds, induced by the presence of a solute molecule. Let 
us recall that, in our model, the presence of too many wa- 
ter molecules close to one another weakens the H bond 
strength, through the c„ parameter, to mimic the fact 
that too low average distance is unfavorable for H bonds. 
On the contrary, a different (lower) weakening parame- 
ter Cs for a solute molecule, can represent a perturbation 
of H bonds, due to an excluded volume effect. Anyway, 
the solute weakening parameter Cs is a way of tuning the 
degree of water ordering induced by the solute. In or- 
der to characterize this effect, together with the role of 
the Cs parameter, let us investigate the H bond coordina- 
tion, that is, the average number of hydrogen bonds per 
molecule. In quite recent papers. Dill and coworkers [l9j 
suggested, on the basis of their off-lattice MB model, that 
this parameter is the appropriate one to rationalize the 
temperature Ts at which the transfer entropy vanishes. 
In particular they distinguished between H bond coor- 
dination for molecules in bulk water and in a hydration 
shell. We can evaluate analogous parameters also for our 
model, in the framework of the first order approximation. 
Each water molecule can form bonds with NN molecules 
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FIG. 3: Site labels for the hexagon probability distribution. 
Thick lines denote triangles for which correlations are taken 
into account. 



only, therefore it is necessary to compute the joint proba- 
bility distributions of a given site with its 6 NNs (hexagon 
probability) . According to the Husimi lattice formulation 
of the first order approximation, it is easy to see that 
only certain triangle correlations are taken into account, 
as shown schematically in Fig|2| Therefore, io being the 
central site configurations and ii . . . the NNs configu- 
rations, the hexagon probability distribution reads 

Pioii...ie — 2 ■ ^ ' 

Pio 

We shall evaluate bond coordinations according to this 
factorization, which, let us notice, is perfectly consistent 
with our approximation scheme [321 ■ Let us consider a 
water molecule in a bonding configuration, for example 
i = 1. It is not necessary to consider also the other 
{i = 2) bonding configuration, due to the fact that the 
liquid phase is isotropic, whereas we do not take into 
account nonbonding configurations (i = 3), because in 
that case the probability of forming a bond is zero. The 
bulk and shell coordinations can be written respectively 
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as 



q=l 
3 

9=1 



(17) 

(18) 



where Pq\Bj or -Pg|5 respectively, is the probability that 
the given molecule forms q bonds, given that its NN sites 
host no solute molecules (bulk water) or at least one so- 
lute molecule (hydration shell). Working in the infinite 
dilution limit, the probability of configurations with more 
than one NN solute is actually a small corrections over 
the probability of having just one solute molecule. Mak- 
ing use of the Bayes theorem, the bulk conditional prob- 
ability can be rewritten as 



q\B 



q,B 



Pe 



(19) 



where Pq^B is the probability that the central molecule 
forms q bonds, and that its NN sites host no solute 
molecules, while Pb is simply the probability that NNs 
host no solute molecules. Moreover, making use also of 
the total probability theorem, the shell conditional prob- 
ability can be rewritten as 



P 



q\S 



q,S 



P, - P 



q.B 



Ps 



1-Pb 



(20) 



where Pq is the total probability that the central molecule 
forms q bonds. 

We can now evaluate the required probabilities, mak- 
ing use of the factorization ((T?)]l . First of all, it is easy to 
see that 



Pq- 



3\ ^12^^1-^12)'"" 



where 



4 4 

Pi = ^ ^ WjWkPljk 
j=0 k=0 
4 

Pl2 = ^WkPl2k- 

k=0 



(21) 

(22) 
(23) 



In fact, the central molecule, in the given i = I configu- 
ration, can form bonds along its 3 arms. Therefore, P12'' 
is the probability that a bond is formed along q given 
arms, while {pi — Pi2)^~'^ is the probability that bond is 
not formed along the remaining 3 — 5 arms. The prob- 
ability factorizes, because arms lie on different triangles, 
which are uncorrelated, according to Eq. Hlt)|) . The bino- 
mial coefficient accounts for different choices of q bonds 
along 3 arms, while the denominator is due to the fact 
that Pq, though not specified by the notation, is a con- 
ditional probability, with respect to the presence of a 
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FIG. 4: H bond coordination for a shell water molecule (solid 
line) and a bulk water molecule (dashed line) vs temperature, 
(a) Cs = 0, (b) Cs = 0.37, and (c) Cs = 0.405. 



central water molecule in the i = I configuration. The 
joint probability Pq^B can be evaluated in similar way 



Pq,B = 



3\ Pu" {Pi - Puf^'' 



pr 



(24) 



where we have to assume that no solute molecule is 
present in the neighborhood, that is 



3 3 

Pi = ^ ^ WjWkPljk 
j=0 k=0 

3 

P12 = ^ WkPi2k = P12 - P124. 
fe=o 



(25) 



(26) 



Finally, the probability that no solute molecule is in the 
neighborhood (bulk condition) can be written as 



Pb = 



Pi 
Pi 



3 • 



(27) 



In the work of Dill and coworkers |l9l| , at low tempera- 
tures H bond coordination for shell water {q)g is greater 
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than H bond coordination for bulk water ((7)3. This be- 
havior is reversed in the high temperature region, where 
H bonds are preferably formed by bulk water. In this 
model, the crossing temperature Tg (where (q) g = (q) ^) 
is approximately equal to Ts (where As* =0). As a 
consequence, average H bond coordination seems to be 
strictly related to the transfer entropy. Let us analyze 
what happens in our model. First of all, wc observe that 
the range of values of H bond coordination is generally 
much lower than the maximum value of 3 H bonds, which 
can be formed by a single water molecule (see Fig. 0J|. 
This fact may be related to a peculiarity of the lattice 
model, in which H bonds can be formed just along given 
directions. Moreover, as far as the temperature behav- 
ior is concerned, upon varying the weakening parameter 
Cs, we can have in principle four different scenarios. As 
shown in Fig.^Ja), at low values of Cs (low H bond weak- 
ening near a solute molecule), shell coordination {q)g is 
greater than bulk coordination {q) ^ for all temperatures. 
Increasing the parameter, it is possible to observe a sce- 
nario similar to the one found by Dill and coworkers: 
a crossing temperature T'g between a low temperature 
region dominated by {q)g and a high temperature one 
dominated by {q)^ [Fig01b)]. Increasing Cs further, a 
more complicated scenario, involving two crossing tem- 
peratures, can be observed. The region where {q)g is 
greater than (g)^ is upper-bounded by Tg, as in the pre- 
vious scenario, and lower-bounded by a different crossing 
temperature [Fig. 2{c)]. This means that the high tem- 
perature behavior is restored at low temperature, though 
in a region very close to the stability limit for the liquid 
phase, which would be very difficult to reach in a real 
system. Finally, for very high values of Cs (high weaken- 
ing of H bonds near a solute molecule) a fourth scenario 
is observed, reversing the behavior of the first scenario, 
that is, {q)g > {q)g for all temperatures. All possible 
behaviors arc summarized in Fig. |S1 where the crossing 
temperatures T5 and T'g arc displayed as functions of the 
weakening parameter Cs- As a comparison, in the same 
figure, we have reported the zero entropy temperature 
Ts, as a function of the same weakening parameter Cj. It 
can be easily observed that, while H bond coordination 
is strongly dependent on the weakening parameter, the 
transfer entropy turns out to be a much more "robust" 
effect, and the T5 temperature has only a relatively slight 
dependence on the weakening parameter. 



V. CONCLUSIONS 

In this paper wc have considered a two-dimensional 
waterlike lattice model, defined on the triangular lat- 
tice, which we had previously shown to reproduce qual- 
itatively typical thermodynamic anomalies of pure wa- 
ter, and extended the model to the case of an aque- 
ous solution. Water molecules are of the Mercedes Benz 
type, with three equivalent bonding arms, while solute 
molecules have no internal degrees of freedom. We have 
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FIG. 5: Weakening parameter Cs vs crossing temperature Tg, 
where {q}g = {q}g, (empty circles interpolated by a solid 
line) and vs Ts, where As* = 0, (filled circles interpolated by 
a dashed line). 



performed our calculations by means of a generalized first 
order approximation on a triangle cluster, which requires 
quite a small computational effort, and had been shown 
to behave quite accurately for the pure water model. 
In particular, in this paper we have addressed the is- 
sue of dilute solutions of inert (apolar) solutes, that is, 
the hydrophobic effect, and we have investigated thermo- 
dynamic equilibrium between liquid and vapor, working 
out solvation quantities in this case. It turns out that the 
model qualitatively reproduces the peculiar features that 
are usually believed to be the fingerprints of hydropho- 
bicity. The solvation free energy is positive, meaning that 
solvation is unfavorable, while the entropy and enthalpy 
display an increasing trend and cross zero at two differ- 
ent temperatures. The solvation heat capacity is large 
and decreases upon increasing temperature. The model 
results compare qualitatively well with experimental re- 
sults about solvation of noble gases into water. We have 
investigated the effect of H bonding, by comparing the 
previously mentioned results with those obtained by set- 
ting to zero the H bond energy parameter. In this case, 
we have obtained transfer quantities that approach the 
ones computed with H bonds at high temperatures, but 
that largely deviates from them upon decreasing temper- 
atures, in the region were H bonding begins to dominate. 
In particular, wc have observed that, while the disaffinity 
of the solute for the solvent remains, and is signaled by 
the positive large value of the solvation free energy, such 
a disaffinity is mainly of enthalpic nature, and both the 
enthalpy and entropy of solvation remain positive at all 
observed temperatures, so that the strong temperature 
dependence, which is typical of the hydrophobic effect, 
disappears. 

In order to a check of the model, we have also investi- 
gated solvation of water into water at liquid- vapor equi- 
librium, for which experimental data are available. We 
have found qualitative agreement in the values of solva- 
tion free energy, entropy and enthalpy, though there is 



9 



some discrepancy in the temperature dependence of en- 
thalphy, which indicates a negative solvation heat capac- 
ity, in contrast with experiments. Though not reporting 
details in the paper, we have verified that this fact is nei- 
ther to be related to the approximation of the ideal/real 
gas phase, nor with the temperature derivative approxi- 
mation, mentioned in Sec. III. On the contrary, we have 
observed that the same kind of discrepancy may be ob- 
served in the case of zero H bond energy, that is, for a 
regular solution. Therefore, we suspect that the discrep- 
ancy is to be related to an intrinsic limitation of lattice 
modelling. The effect is relatively small, so that it is com- 
pletely invisible, when the dominant effect of H bonding 
is turned on. 

Finally, we have performed a calculation of the aver- 
age number of H bonds formed by a single water molecule 
(H bond coordination, in short), in two different cases, 
namely when the molecule is placed in the first hydra- 
tion shell of a solute molecule (shell coordination), and 
when it is not (bulk coordination). According to Dill's 
Mercedes Benz model, these two parameters are clearly 
related to the solvation entropy behavior, namely, neg- 
ative solvation entropies (low temperature) correspond 
to higher shell coordination, while positive solvation en- 
tropies correspond to higher bulk coordination. Due to 
the fact that our purpose was to obtain a kind of lattice 



version of the MB model, describing essentially the same 
physics underlying the observed phenomenology, we have 
tried to verify whether the same effect could be observed 
in our model or not. The answer is essentially no, but 
some interesting observations can be done. We have per- 
formed the investigation, upon varying the solute weak- 
ening parameter, which, in our model, is a way of tuning 
the degree of water ordering induced by the solute. We 
have observed that such parameter strongly affects the 
behavior of H bond coordination, and in particular the 
temperature at which shell and bulk coordinations be- 
come equal. On the contrary the temperature of zero en- 
tropy Ts, which is actually one of the striking features of 
the hydrophobic effect, is quite "robust" and relatively 
independent of variations of the weakening parameter. 
Two questions arise from the discussed behavior. On the 
one hand, we may suspect that the lattice model is def- 
initely too simple to capture the microscopic physics of 
the hydrophobic effect in a correct way, or that the ap- 
proximation level is insufficient. But the robustness of 
the zero entropy temperature may also suggest that the 
simple relationship between the balance of bulk and shell 
H bonds and the zero of transfer entropy, proposed by 
Dill and coworkers, may be model-dependent, and ought 
to be verified more carefully. We plan to investigate on 
such issues in future works. 
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